{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# A Cantera Simulation with Reaction Sensitivity: Comparison with Native RMG Sensitivity Analysis and CHEMKIN Sensitivity Analysis\n",
    "\n",
    "\n",
    "Note that this requires Cantera with SUNDIALS installed for sensitivity support.  If you are using Anaconda,\n",
    "Cantera version >= 2.3.0 is required"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's check what version of Cantera you have installed"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import cantera\n",
    "print cantera.__version__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from rmgpy.chemkin import *\n",
    "from rmgpy.tools.canteraModel import *\n",
    "from rmgpy.tools.plot import parseCSVData\n",
    "from rmgpy.species import Species\n",
    "from IPython.display import display, Image"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load the species and reaction from the RMG-generated chemkin file `chem_annotated.inp` and `species_dictionary.txt` file found in your `chemkin` folder after running a job."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "speciesList, reactionList = loadChemkinFile('data/minimal_model/chem_annotated.inp',\n",
    "                                            'data/minimal_model/species_dictionary.txt',\n",
    "                                           'data/minimal_model/tran.dat')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Set the reaction conditions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Find the species: ethane and methane\n",
    "user_ethane = Species().fromSMILES('CC')\n",
    "user_methane = Species().fromSMILES('C')\n",
    "speciesDict = getRMGSpeciesFromUserSpecies([user_ethane, user_methane], speciesList)\n",
    "ethane = speciesDict[user_ethane]\n",
    "methane = speciesDict[user_methane]\n",
    "sensitiveSpecies = [ethane, methane]\n",
    "\n",
    "#reactorTypeList = ['IdealGasReactor']\n",
    "reactorTypeList = ['IdealGasConstPressureTemperatureReactor']\n",
    "molFracList=[{ethane: 1}]\n",
    "Tlist = ([1300],'K')#,1500,2000],'K')\n",
    "Plist = ([1],'atm')\n",
    "reactionTimeList = ([0.5], 'ms')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Create cantera object, loading in the species and reactions\n",
    "job = Cantera(speciesList=speciesList, reactionList=reactionList, outputDirectory='temp', sensitiveSpecies=sensitiveSpecies)\n",
    "# The cantera file must be created from an associated chemkin file\n",
    "\n",
    "# We can either load the Model from the initialized set of rmg species and reactions\n",
    "job.loadModel()\n",
    "\n",
    "# Or load it from a chemkin file by uncommenting the following line:\n",
    "#job.loadChemkinModel('data/minimal_model/chem_annotated.inp',transportFile='data/minimal_model/tran.dat')\n",
    "\n",
    "# Generate the conditions based on the settings we declared earlier\n",
    "job.generateConditions(reactorTypeList, reactionTimeList, molFracList, Tlist, Plist)\n",
    "# Simulate and plot\n",
    "alldata = job.simulate()\n",
    "job.plot(alldata)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Show the plots in the ipython notebook\n",
    "for i, condition in enumerate(job.conditions):\n",
    "    print 'Cantera Simulation: Condition {0} Mole Fractions'.format(i+1)\n",
    "    display(Image(filename=\"temp/{0}_mole_fractions.png\".format(i+1)))\n",
    "    \n",
    "    print 'Cantera Simulation: Condition {0} Ethane Reaction Sensitivity'.format(i+1)\n",
    "    display(Image(filename=\"temp/{0}_ethane(1)_sensitivity.png\".format(i+1)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Let's compare against the same simulation in RMG\n",
    "# Create an input file\n",
    "\n",
    "input = '''\n",
    "database(\n",
    "    thermoLibraries = ['primaryThermoLibrary'],\n",
    "    reactionLibraries = [],\n",
    "    seedMechanisms = [],\n",
    "    kineticsDepositories = 'default',\n",
    "    kineticsFamilies = 'default',\n",
    "    kineticsEstimator = 'rate rules',\n",
    ")\n",
    "\n",
    "species(\n",
    "    label = \"ethane\",\n",
    "    structure = SMILES(\"CC\"))\n",
    "\n",
    "species(\n",
    "    label = \"methane\",\n",
    "    structure = SMILES(\"C\"))\n",
    "\n",
    "simpleReactor(\n",
    "    temperature = (1300,\"K\"), \n",
    "    pressure = (1,\"atm\"),\n",
    "    initialMoleFractions={\n",
    "        \"ethane\": 1\n",
    "    },\n",
    "    terminationTime = (0.5,\"ms\"),\n",
    "    sensitivity=['ethane','methane']\n",
    ")\n",
    "\n",
    "model(\n",
    "    toleranceMoveToCore = 0.04,\n",
    ")\n",
    "\n",
    "options(\n",
    "    saveSimulationProfiles = True,\n",
    ")\n",
    "\n",
    "'''\n",
    "f = open('temp/temp_input.py', 'wb')\n",
    "f.write(input)\n",
    "f.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from rmgpy.tools.sensitivity import runSensitivity\n",
    "runSensitivity('temp/temp_input.py', 'data/minimal_model/chem_annotated.inp', 'data/minimal_model/species_dictionary.txt')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print 'RMG Native Simulation: Species Mole Fractions'\n",
    "display(Image(filename=\"temp/solver/simulation_1_19.png\"))\n",
    "\n",
    "print 'RMG Native Simulation: Ethane Reaction Sensitivity'\n",
    "display(Image(filename=\"temp/solver/sensitivity_1_SPC_1_reactions.png\"))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Let's also compare against the same simulation and sensitivity analysis that was conducted in CHEMKIN\n",
    "# and saved as a .csv file\n",
    "time, dataList = parseCSVData('data/minimal_model/chemkin_mole_fractions.csv')\n",
    "SimulationPlot(xVar=time,yVar=dataList).plot('temp/chemkin_mole_fractions.png')\n",
    "print 'CHEMKIN Simulation: Mole Fractions'\n",
    "display(Image(filename=\"temp/chemkin_mole_fractions.png\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "time, dataList = parseCSVData('data/minimal_model/chemkin_sensitivity_ethane.csv')\n",
    "ReactionSensitivityPlot(xVar=time,yVar=dataList).barplot('temp/chemkin_ethane_sensitivity.png')\n",
    "print 'CHEMKIN Simulation: Ethane Reaction Sensitivity'\n",
    "display(Image(filename=\"temp/chemkin_ethane_sensitivity.png\"))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
